Landau functions for non-interacting bosons 
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We discuss the statistics of Bose- Einstein condensation (BEC) in a canonical ensemble of A'^ non- 
interacting bosons in terms of a Landau function C^'~'{q) defined by the logarithm of the probability 
distribution of the order parameter q for BEC. We also discuss the corresponding Landau function 
for spontaneous symmetry breaking (SSB), which for finite A'^ should be distinguished from C^'^{q). 
Only for N ~> oo BEC and SSB can be described by the same Landau function which depends on 
the dimensionality and on the form of the external potential in a surprisingly complex manner. 
For bosons confined by a three-dimensional harmonic trap the Landau function exhibits the usual 
behavior expected for continuous phase transitions. 

PACS numbers: 03.75.Hh, 05.70.Fh 



The concept of a Landau function C{p) has been ex- 
tremely useful to understand the nature of phase transi- 
tions 1]. For example, in the vicinity of a second order 
(continuous) thermal phase transition C{p) smoothly de- 
velops minima at non-zero values of the relevant order 
parameter p when the temperature T drops below the 
critical temperature T^, see Fig. Q] (a). On the other 
hand, in the case of a first order phase transition the lo- 
cation of the minima of C{p) changes discontinuously at 
T = Tc, as shown in Fig. ^b). While for infinite systems 
and in the absence of spontaneous symmetry breaking 
(SSB) C{p) can be identified with the Gibbs potential, 
this identification is in general not valid for finite sys- 
tems or in the symmetry broken phase. In the latter case 
the Gibbs potential, which by construction is a convex 
function of p, exhibits a plateau as shown in Fig. ^ Mi- 
croscopically, C{p) is defined in terms of the logarithm 
of the probability distribution of the order parameter in 
a given statistical ensemble Q. From this probabilistic 
definition it is clear that C{p) does not have a strictly 
thermodynamic interpretation. 

For non-interacting systems the Landau function is in 
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FIG. 1: Typical evolution of the Landau function £ as a func- 
tion of the relevant order parameter p close to (a) second order 
and (b) first order phase transitions. The dashed line denotes 
the plateau of the corresponding Gibbs potential for T <Tc. 



most cases not very interesting, because the correlations 
responsible for phase transitions are usually due to inter- 
actions between the microscopic degrees of freedom. An 
important exception are non-interacting bosons, where 
at low temperatures the correlations imposed by quan- 
tum statistics are sufficient to give rise to a phase tran- 
sition, the Bose- Einstein condensation (BEC) Q,!^. In 
the condensed state the fraction of particles that occupy 
the single-particle state with the lowest energy is of the 
order of unity. More precisely, the order parameter for 
BEG can be chosen to be the expectation value of the 
operator q — bj^bo/N, where 6o annihilates a boson in the 
single-particle ground state and N is the total number of 
particles in the system. 

Conceptually, BEG should be distinguished from SSB, 
which manifests itself via a finite expectation value of the 
operator (f) = b^/^fN in the limit iV oo. Of course, 
for any finite TV there is no off-diagonal long-range or- 
der, so that (0) = 0. For interacting bosons in the limit 
^ oo it can be shown that {q) — whenever 
the Bogoliubov approximation of replacing the operator 
6o by a complex number yields the right pressure 
For finite systems, however, one should distinguish be- 
tween SSB and BEG. More precisely, both phenomena 
should be characterized by different probability distribu- 
tions, Pf^^[(j)) Etnd P^^^(g), which depend on the eigen- 
values q and 4> of the operators q and 4> defined above. 
We parameterize the probability distributions in terms of 
Landau functions C^^^{(f)) and £^^^(g), see Eqs. ((71g|) 
below. For interacting bosons the condensation transi- 
tion is known to be second order 0, so that the Landau 
functions resemble in this case Fig. ^(a). But how do 
the Landau functions of non-interacting bosons look like? 
Because the order parameter is continuous, the first or- 
der scenario in Fig.n(b) can be ruled out. In spite of the 
enormous activity in the field of BEC in recent years, this 
question has apparently not been addressed in the litera- 
ture. Although the probability distribution P^-^ p (q ) has 
been studied by several authors [3, IE la la lid, the 
associated Landau function £^^*-^(g) has not been con- 
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sidered. The corresponding quantities for SSB have not 
even been properly defined in the hterature. 

To describe BEC experiments in finite clusters of 
bosonic atoms [Tl| . we use a statistical ensemble with 
constant particle number TV, i.e., the canonical en- 
semble. The mathematically more convenient grand 
canonical ensemble is neither relevant for experiments 
\n\ . nor is it appropriate to describe the condensed 
phase 0,0 • To derive explicit expressions for £^^(0) 
and C^^{q), consider the canonical partition function, 
Z]\j = Tr^v exp[— /37J], where (3 = 1/T is the inverse tem- 
perature and Trjv denotes the trace over the A^-particle 
Hilbert space of a system of non-interacting bosons with 
Hamiltonian H ~ J2m ^■mblnbm- Here, bm annihilates 
a boson in a single particle state with quantum num- 
bers m and energy Em- For free bosons (FB) with 
mass M confined to a _D-dimensional volume with 
periodic boundary conditions Em = (27r?im/L)^/2M, 
where m is a Z?-dimensional vector with integer compo- 
nents rrii — 0, ±1, ±2, . . ., i — 1, . . . , D. For bosons con- 
fined by a harmonic potential (HB) with oscillation fre- 
quency Lo the energy is Em — fiuj + 1/2), where 
nii = 0, 1, 2, . . .. For later convenience we introduce di- 
mensionless energies em — PEm and the relevant dimen- 
sionless density p = {Xth/ L)^ N , where for FB the ther- 
mal de Broglie wavelength is Ath = ?iy27r/3/M, whereas 
for HB we define Ath = huj^L and L = {Muj/hy/^. 

With the help of the projection operator 5^ ^ — 

Jo" l^e'y(''-^\ where N^J^^ bl^bm, we can express 
Zn in terms of a trace over the entire Fock space 



Zn = Tr 
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The partial trace over the m 7^ sector of the Fock space 
yields 

Tr^^oe"^'"^''<''"+*^)^^'''" = e-^^'^'^^'y''^^ , (2) 

where 7 = p/N and Joia,!) = lY.m^o^'^l'^ ~ e""'"^"]- 
The contribution from the m = sector can be written 
in terms of coherent states \z) satisfying ^ok) = ^k) a-s 
follows, 
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where cPz = rfRezrflmz, and in the second line we have 
evaluated the matrix element by inserting the resolution 
of unity in the basis \n) of particle number eigenstates, 
b\bQ\n) — n\n). With (j) — zj\fN we finally obtain 



Zn = d (pe 



(4) 
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FIG. 2: (Color online) Order parameter probability distribu- 
tions (left column) and Landau functions (right column) of 
= 10 bosons in a three-dimensional harmonic potential for 
difTerent dimensionless densities p — {Xth/ L)'^ N . The plots in 
the second row correspond to the critical density pc « 2.85 
for = 10. Discrete points denote P^^'"' (plotted versus 
q = n/N) and C^'~' (versus y^), while the continuous curves 
represent P^^^ (versus 10^) and Cjf^ (versus Re^). 



where the Landau function for SSB is 



N 



(0) = i0i -^H^l^^^ 



and the corresponding Landau function for BEC is 



(5) 
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We have normalized the Landau functions such that they 
are dimensionless and have finite limits for iV — > 00. We 
conclude that the normalized probability density of the 
order parameter (p for SSB is 



Pf^(0) = Z^ie-^^«"W), 



(7) 



while the probability of observing a fraction q = n/N of 
bosons in the single-particle state with lowest energy is 



3BEC 

N 



(8) 



For HB the probability P^^'~^{q) has been calculated nu- 
merically for N < 500 by Balazs and Bergeman ^ using 
an efficient recursive algorithm; however, they neither 
considered the corresponding Landau function nor the 
analogous quantities for SSB. 

In Fig. 121 we show numerical results for the above Lan- 
dau functions and the associated probability distribu- 
tions for HB in D = 3. Obviously, the behavior is typical 
for second order phase transitions. We shall show shortly 
that for HB in 13 = 3 this remains true for ^ 00. In 
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this limit the critical (dimensionless) density for BEG 
and SSB in D dimensions is Pc = Cid), where d — D/2 
for FB and d = Z? for HB. For finite clusters with N > 2 
particles we define the critical density pf^^{N) for BEG 
as the value of p where ^^.'^^(O) = C^^'=^{l/N). Analo- 
gously, at p = Pc^^iN) the minimum of £|p^(0) shifts 
from = to 101 > 0. This definition of the critical den- 
sity (and the associated critical temperature) in finite 
clusters is free of ambi guit ies and remains valid for inter- 
acting bosons IEIEBEI- Of course, in finite systems 
one may use other criteria, such as the behavior of the 
specific heat 0], to define a characteristic temperature 
which approaches the critical temperature for N —^ oo. 
Numerically we find that even for small clusters contain- 
ing a few bosons pf'^'^{N) « p^^^{N). The most efficient 
way to obtain the Landau functions is to directly calcu- 
late P^^'^{q) recursively 0, and then use Eqs. (|5I7I8|I . 
As expected, with increasing N the critical density slowly 
approaches the bulk value pc = C(3) ~ 1-09 from above 
(for example pc(lOO) « 1.93 and pc(500) « 1.60), and the 
probability distributions develop narrows peaks, so that 
the difference between the most probable and the average 
value of the order parameter vanishes. The advantage of 
parameterizing the probability distributions in terms of 
Landau functions is that only their detailed shape (but 
not their overall scale) depends on N, because the lead- 
ing A^-dependence has already been taken into account 
via the prefactor N in the exponent exp[— A>Cjv]- 

In the rest of this work we consider the limit N oo 
with fixed p. For FB this is the usual thermodynamic 
limit, while for HB this limiting procedure amounts to 
letting — > at fixed uj-'^N oc p, see Ref. Then the 
summation in Eq. ((SJ is sharply peaked at n = A^ji^p, so 
that £L'i^('/')_= ^^''m') ^ C^im where C^{q) = 
— limjv_,oo N ^ ln/jv(<z), and /Ar(g) is the following inte- 
gral in the complex z ^ x + iy plane. 



(9) 



The complex function f(z) is in D dimensions given by 

f{z) = zil-q) + p-'gd+i{z) . (10) 

RecaU that d ^ D/2 for FB and d = D for HB. The 
Bose-Einstein integral gs{z) is defined by 
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du- 



'J.+Z _ 2 



(11) 



were T(s) is the F-function. In the complex z-plane gs{z) 
has infinitely many branch points at integer multiples of 
27rj, with branch cuts parallel to the negative real axis. 

The leading asymptotic behavior of the In {q) for large 
N can be obtained using the method of steepest descends 
[ill. For p < p^= ga{0) = C(rf), or if p > pc but q > 
qc = 1 — Pc/ P > 0, the saddle point equation 
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FIG. 3: (Color online) Typical behavior of'Ref{x+iy) defined 
in Eq. 1)10^ in the regime q > q^. The graph is for d = 3/2, 
q = 0.95, and qc = 0.05. The end points of the integration 
and the saddle points are marked by black dots. The thick 
line marks the value of Re/ on a deformation of the original 
integration contour in Eq. (|5J along steepest descend paths 
through the saddle points zq — a{q, p) and zi = a{q, p) + 2ni. 



has a solution zq — a{q, p) on the positive real axis. Using 
the analyticity of the integrand, we may then deform the 
original integration contour into a combination of steep- 
est descend paths through zq and the related saddle point 
zi = a{q, p) + 27ri, as shown in Fig. |31 The asymptotic 
behavior of 1^ (q) for A^ — > oo is determined by the value 
of the integrand at the saddle points zq and zi, which 
yield complex conjugate contributions. We thus obtain 
in this regime Loo{q) — — f {oi{q, p)) . The correspond- 
ing effective potential C/eff(0) — -Cood^P) — 'Coo(O) shows 
rather non-trivial behavior, in spite of the fact that we 
are dealing with a non-interacting system. For simplicity, 
we consider here only the regime close to the phase tran- 
sition, where \p — Pc\ Pc- In the normal phase p < pc 
we may expand UcS in powers of 
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Here ao = a{q — 0, p) is the real solution of Eq. H12|l 
for q = 0, i.e. p = gd{ao), and uq — X^^('^o)) where 
x(q:) = p^^gd-i{a) is the susceptibility of the Bose gas. 

Eq. (|13|l is the usual behavior of the effective potential 
in the vicinity of continuous phase transitions. For p 
Pc the parameter ao vanishes, so that at the first sight 
it seems that the leading term of the effective potential 
at the critical point is proportional to |(/)|'*. However, 
this conventional scenario is only correct for d > 2, i.e. 
for dimensions D > i (FB) or D > 2 (HB). For d < 
2 the susceptibility xi^t) diverges for a ^ 0, implying 
uo — > for p ^ Pc- For 1 < d < 2 we find at the 
critical point to leading order C/off(0) ~ A{Pc)\4>\'^ with 
A{p) = (l-d"i)|p/r(l-d)|3^. In particular, for FB in 

D = 3 we obtain J7eff((/)) = In the condensed 

phase p > pc the Landau function for q > q^ is still 
determined by the above saddle points zq and zi . In this 
case we obtain close to the phase transition 



df{z)/dz = 1 - gr - p ^gd{z) = 
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for 2 > d > 1 
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FIG. 4; (Color online) Kef{x + iy) and steepest descend paths 
for q < Qc- The graph is for d — 3/2, Qc = 0.05, and q — 0.045. 

Note that for q > Qc the Landau function can be iden- 
tified with the Gibbs potential for ^ oo, so that the 
above results can also be derived from thermodynamics. 

On the other hand, in the condensed phase p > pc the 
Gibbs potential is constant for q < qc — I — Pc/ P, and 
is not necessarily equal to the Landau function Ij. It 
turns out that in this regime the saddle point equation 
(|12|l does not have solutions on the principal sheet of the 
Riemann surface of the complex function f{z) defined 
in Eq. However, there are saddle points on higher 

sheets, two of which (which we call zq and zi) are the ana- 
lytic continuation of the saddle points zq and zi discussed 
above for q > qc- While we have not attempted a global 
analysis of the descend paths on the higher sheets of the 
Riemann surface of f{z) for arbitrary d, we can show 
that for d > 5/3 the path connecting the origin to zq 
continues to be a steepest descend path, so that the sad- 
dle points Zq and zi are still relevant for the asymptotics 
of In{(1) for — > oo. The leading behavior of Coo{q) for 
small qc ~ q > can then be obtained by analytic con- 



tinuation of the corresponding result for q > q^ which 
amounts to replacing |<^p — qc by qc — \4>\'^ in Eq. H14|l . 
Hence, for d > 5/3 the effective potential qualitatively 
resembles the "Mexican hat" characteristic for continu- 
ous phase transitions, although for 5/3 < d < 2 its shape 
cannot be approximated by a quartic potential |17| |. 

On the other hand, for 1 < d < 5/3 the saddle points 
Zq and zi become inadmissible for the asymptotic ex- 
pansion of lN{q) in the sense defined on p. 267 of Ref. 
[iq : the end-points of the integration can then be di- 
rectly connected by steepest descend paths avoiding any 
saddle point. For d = 3/2 the proper deformation of 
the integration contour is shown in Fig. 0] As a con- 
sequence, for d < 5/3 the Landau function is constant 
Coo{q) = -p~^gd+iiO) for < g < (?c = 1 - Pc/p, just 
like the Gibbs potential. Note that in I? = 3 the behavior 
of Cryo{q) depends on the form of the external potential: 
while for FB {d — D/2 — 3/2) the Landau function has 
a plateau for < g < gc, for HB {d — D — 3) it has the 
"Mexican hat" form. 

In summary, we have defined and evaluated the Lan- 
dau functions and the order parameter probability distri- 
butions for BEG and SSB of non-interacting bosons in a 
canonical ensemble. The shape of the Landau functions 
in the condensed phase depends on the dimensionality 
and on the shape of an external potential in a surprisingly 
complex way. For bosons in a harmonic trap the evolu- 
tion of the Landau functions with density or temperature 
is typical for second order phase transitions. The shape 
of Landau functions leads to the most natural extension 
of the concept of the critical temperature for BEG and 
SSB in finite systems. 
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